Two and three-dimensional oscillons in nonlinear Faraday resonance 
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We study 2D and 3D localised oscillating patterns in a simple model system exhibiting nonlinear 
Faraday resonance. The corresponding amplitude equation is shown to have exact soliton solutions 
which are found to be always unstable in 3D. On the contrary, the 2D solitons are shown to be stable 
in a certain parameter range; hence the damping and parametric driving are capable of suppressing 
the nonlinear blowup and dispersive decay of solitons in two dimensions. The negative feedback loop 
occurs via the enslaving of the soliton's phase, coupled to the driver, to its amplitude and width. 



o 
o 

(N 
o 

Q 

oo 



00 
Oh 

d 



> 
o 

(N 



O 

15 

> 

X 



Oscillons are localised two-dimensional oscillating 
structures which have recently been detected in exper- 
iments on vertically vibrated layers of granular mate- 
rial 0, Newtonian fluids and suspensions Numer- 
ical simulations established the existence of stable oscil- 
lons in a variety of pattern- forming systems, including 
the Swift-Hohenberg and Ginsburg-Landau equations, 
period-doubling maps with continuous spatial coupling, 
semicontinuum theories and hydrodynamic models [^|,[| . 
Although these simulations provided a great deal of in- 
sight into the phenomenology of the oscillons (in particu- 
lar, demarcated their existence area on the corresponding 
phase diagrams), little is known about the mechanism by 
which they acquire or loose their stability. 

In this Letter, we consider a model equation which has 
exact oscillon solutions and allows an accurate charac- 
terisation of their existence and stability domains. The 
main purpose of this work is to understand how the 
oscillons manage to resist the general tendencies to- 
ward nonlincarity-induced blow-up or dispersive decay 
which are characteristic for localised excitations in two- 
dimensional media. Our model admits a straightforward 
generalisation to three dimensions and we use this op- 
portunity to explore the existence of stable oscillons in 
3D as well. 

The model consists of a Z?-dimensional lattice of para- 
metrically driven nonlinear oscillators (e.g. pendula) || 
with the nearest-neighbour coupling: 
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Assuming that the coupling is strong: k ^> 1; that the 
damping and driving are weak: a — 7s 2 , p = 2he 2 where 
e <C 1; and that the driving half-frequency is just below 
the edge of the linear spectrum gap: uj 2 = 1 — e 2 , the os- 
cillators execute small-amplitude librations of the form 
k = 2etp(t,x k )e- 

Xk 



c.c. + 0(e 3 ), where t = £ 2 t/2, 



the parametrically driven damped nonlinear Schrodinger 
(NLS) equation. In 2D, this equation was invoked as a 
phenomenological model of nonlinear Faraday resonance 
in water j|. It also describes an optical resonator with 
different losses for the two polarisation components of 
the field ||. In the absence of the damping and driv- 
ing, all localised initial conditions in the 2D and 3D NLS 
equation are known to either disperse or blow-up in fi- 
nite time ^) . Surprisingly, numerical simulations of (Q) 
with sufficiently large h and 7 revealed the occurrence of 
stable (or possibly long-lived) stationary localised excita- 
tions ||]. However no analytic solutions were found, and 
a possible stabilisation mechanism remained unclear. 

In fact there are two exact (though not explicit) sta- 
tionary radially-symmetric solutions given by 



A ± e- W ±n Q (A ± r); (r 2 



(3) 



1 ± ^Jh 2 - 7 2 , 9+ = iarcsin( 7 //i), 6L 



where A\ 

% — 6 + , and 7£o(r) is the bell-shaped nodeless solution of 

V 2 r K - K + 2V? = 0; TZ r (0) = ft (00) = 0. (4) 

(Below we simply write ft for ft .) In (Q), V 2 = 
d 2 + (D - l)r- x a r . Solutions of Eq.(§) in D = 2 and 
3 are well documented in literature. (See e.g. and 
refs therein.) One advantage of having an explicit de- 
pendence on h and 7, is that the existence domain is 
characterised by an explicit formula. The soliton ex- 
ists for all h > 7; the ip~ exists for 7 < h < y/l + -f 2 . It 
is pertinent to add here that for h < 7, all initial condi- 
tions decay to zero. This follows from the rate equation 

dtM 2 = l[r(XrM 2 )r}r + 2\^\ 2 (h S m2 X ~ 7), (5) 

where if) = \ip\e~ lx . Defining N — J \ip\ 2 dx, Eq.(|^) im- 
plies N t < 2(h - j)N whence N(t) -> as t -> 00. 

We now examine the stability of the two solitons. Lin- 
earising Eq.(|^) in the small perturbation 



k and the slowly varying amplitude satisfies 



<ty>(x, t) = e (^-r)t-^± [ u (x) + i«(x)], 



(6) 



itpt + VV -I- 2\^\ 2 ip -ip = hip* - 171/", 



(2) where x = A±x, t = A\t, we get an eigenvalue problem 
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Liu = -(ji + T)v, (La - e)v = (ji - F)u, (7) 

where T — j/A± and the operators 

L = -V 2 + l-2K 2 (r), L x = L - 4K 2 (r), (8) 

with V 2 = Yli=i d 2 1 ■ (We are dropping tildas be- 
low.) The quantity e, e = ±2^/h 2 — j 2 /A±, is positive 
for the ip + soliton and negative for . Each e defines a 
"parabola" on the (h, 7)-plane: 



h= v/e 2 /(2-e) 2 + 7 2 



(9) 



Introducing A 2 — /i 2 — T 2 and changing v(x) — > (/i + 
T)A- 1 w(x) (TJ, Eq.(@) is reduced to a one-parameter 
eigenvalue problem: 



(Lo — e)i> = Alt, Liu 



-Am. 



(10) 



Since lZ (r) is nodeless in < r < oo, and io^-o = 0, 
the operator L — e is positive definite for e < 0. In this 
case the eigenvalue can be found as a minimum of the 
Rayleigh quotient: 



- A = min • 



(w\L\\w) 
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w 



(11) 



The operator L\ has D zero eigenvalues associated with 
the translation cigcnfunctions d{R.(r), i = 1,2, ...D; 
hence it also has a negative eigenvalue with a radial- 
symmetric eigenfunction wa(r). Substituting wq into the 
quotient in (|ll|), we get — A 2 < whence fi > T. Thus the 
soliton %p~ is unstable (against a nonoscillatory mode) for 
all D, h and 7, and may be safely disregarded. 

Before proceeding to the stability of (for which we 
have e > 0), we make a remark on the undamped, un- 
driven case (e = 0.) In 3D, the eigenvalue problem (M) 
has a zero eigenvalue associated with the phase invari- 
ance of the unperturbed NLS equation (Q) and another 
one, associated with the scaling symmetry: 



L 
Li 



n 

-\(rH) 7 



(12) 



Both the eigenvector (1Z, 0) T and the rank-2 generalised 
eigenvector (0, — ^(rlZ) r ) T are radially-symmetric. In 
2D the number of repeated zero eigenvalues associated 
with radially-symmetric invariances is four; in addition 
to those in (jl^) we have a two-parameter group of the 
lens transformations giving rise to 
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(13) 



with some g(r). When h 2 — 7 2 (or, equivalcntly, e) devi- 
ates from zero, all the above invariances break down and 
the two (respectively, four) eigenvalues move away from 
the origin on the plane of complex A. The directions of 
their motion are crucial for the stability properties. 



We can calculate A(e) perturbatively, assuming 

A = Aie^ + \ 2 ei + A 3 e^ + (14) 



u = uie 4 + 1A26 4 



v = 7Z + vie 4 + V2t l + 



where = i>j(r), Ui — u.i(r). Substituting into (|l(i|), the 
order e 1 / 4 gives U\ = —\iL^ l 1Z. Using (|l2|), U\ is found 
explicitly: u\ = (\i/2)(rlZ) r . At the order e 2 / 4 we get 
U2 = —^L^IZ. and equation L0V2 = XiUi. Since Lq has 
a null eigenvector, Tt(r), this equation is only solvable if 



Ai / K{r) Ul (r)d-x 



-^l^f 1 J ^V)dx = 0. (15) 



In the two-dimensional case the condition < \15\ ) is satisfied 
for any Ai whereas in D = 3 we have to set Ai = 0. Next, 
at the orders e 3 / 4 and e 4 / 4 we obtain, respectively, 



L v 3 = X 2 ui + Xiu 2 = XiX 2 (rTZ) r , 

L V4 = R + A1M3 + A2M2 + A3W4. 



(16) 
(17) 



Eq.(|16|) is solvable both in 2D and 3D. The solvability 
condition for ( |l7|) reduces to 

fn 2 d x 



At 



(TZ\L^L^L^\n) 



A 



16 jTZ 2 r 2 d^ 



{K\L^\1l) 



(18) 
(19) 



in two and three dimensions, respectively. 

Thus we arrive at two different bifurcation scenarios. 
In 3D, where Ai = and A2 is real, two imaginary eigen- 
values ± I A2 1 e 1 / 2 converge at the origin as e — > from the 
left. (This does not mean that the -0 _ soliton is stable as 
there still is a pair of finite real eigenvalues for e < 0.) As 
e grows to positive values, the imaginary pair ±|A 2 |e 1 / 2 
moves onto the real axis. A numerical study jll|] of the 
eigenvalue problem ( |l0| ) shows that when e is further in- 
creased, the four real eigenvalues collide, pair wise, and 
acquire imaginary parts. Importantly, for all < e < 1 
the imaginary parts remain smaller in magnitude than 
the real parts. This means that Re^t remains greater 
than r all the time, implying that the three-dimensional 
4> + soliton is unstable for all h and 7. 

The bifurcation occurring in 2D is more unusual. As e 
approaches zero from the left, four eigenvalues converge 
at the origin, two along the real and two along imagi- 
nary axis: A sa ±|Ai|(— e) 1 / 4 , ±i|Ai|(— e) 1 / 4 . As e moves 
to positive, the four eigenvalues start diverging at 45° to 
the real and imaginary axes. Hence to the leading order, 
ImA s» Re A, and in order to make a conclusion about 
the stability, we need to calculate the higher-order cor- 
rections. The order e 5 / 4 produces a solvability condition 



xiX2(n\L^ L L^ L L^ L \n) 



a|A2 

16 



K 2 r 2 d* = 0, 



yielding A2 = 0. (Here we made use of (|l3|).) Finally, the 
order e 6 / 4 defines A3 (where g(r) is as in (|13|)): 
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FIG. 1. Stability diagram for two-dimensional solitons. The 
(7, h — 7)-plane is used for visual clarity. No localised or periodic 
attractors exist for h < 7 (below the horisontal axis). The region 
of stability of the soliton ip + lies to the right of the solid curve. 
The dotted curve gives the variational approximation to the sta- 
bility boundary of the %p + soliton: h = (1 + 7 4 "! 1 / 2 



A, = -L 



A? / g(r)K(r)r 2 dyL 



Ai 2 jTZ 2 (r)r 2 dx 



(20) 



Taking Ai in the first quadrant, Ai = e l7r / 4 |Ai|, and 
doing the integrals in ( |l8| ) , (|2(]) numerically, we conclude 
that A3 is in the second quadrant, A3 = e 3l7T / 4 |A3|, which 
implies that |ImA| > |ReA|. In terms of A, the stability 
criterion Re^t < T is written as 7 > j c , where 



7c(e) 



ReA(e)ImA(e) 



2-e ^(ImA) 2 - (ReA) 2 



(21) 



The smallest 7 for which the soliton can be stable, is 
given by 



lim 7c ( e ) = -i=|A 1 | 3 / 2 |A3r 1 / 2 . 
2V2 



(22) 



Substituting for Ai, A3 their numerical values, ( p2| ) gives 
7c(0) = 1.00647. For e / we obtained A(e) by solving 
the eigenvalue problem (10) directly Here we have 
restricted ourselves to radially- symmetric u(r) and v(r). 
Expressing e via 7 C from (f2l| ) and feeding into (||) , we get 
the stability boundary on the (h, 7)-plane (Fig.l). 

Asymmetric perturbations do not lead to any instabil- 
ities in 2D. To show this, we factorise, in (|l~0|), u(x) = 
u(r)e zmtp and v(x) = v(r)e lmtp , where tantp = y/x and m 
is an integer. The eigenproblem (|h]) remains the same, 
with only the operators Lq and L\ being replaced by 



4 m) EEL + m 2 /r 2 , 



• (m) 



Li + m 2 /r 2 



(23) 



The crucial observation now is that L ( m > with m 2 > 1 
does not have any (not even positive) discrete eigenval- 
ues. Wc verified this numerically for m 2 = 1; this rules 



out their appearance for all other to. Therefore the op- 
erator £q — e with e < 1 is positive definite, and the 
eigenvalues of the problem (|ic| ) can be found from the 
variational principle (0). The operator has a zero 
eigenvalue with the eigenfunction w^'(r) = TZ r (f) which 
has no nodes for < r < 00; hence its all other eigen- 
values (if exist) are positive. This also implies that L^ 1 
with m 2 > 1 are positive definite. Thus the minimum of 

1 and positive 



the Rayleigh quotient (|11| ) is zero for m 2 
for m 2 > 1. 

Besides the nodeless solution lZo(r), the "master" 
equation (^) has solutions 1Z n (r) with n nodes, n = 
1,2,.... These give rise to a sequence of nodal solu- 
tions of the damped-driven NLS (||), defined by Eq.(||) 
with IZq — * lZ n . It is easy to realise that the solitons 
ip~ are unstable against radially-symmetric nonoscilla- 
tory modes for all h, 7, n and D. (The proof is a simple 
generalisation of the one for ipQ .) To examine the sta- 
bility of the ■f/'n soliton, we solved the eigenvalue prob- 
lem ([n]) numerically, with the operators L ™^ as in (|23|). 
In 3D, positive real eigenvalues (with radially-symmetric 
eigenfunctions) are present in the spectrum for all e; thus 
the three-dimensional nodal solitons are always prone to 
a symmetric collapse or dispersive spreading. In 2D, the 
ipn solitons are stable against radially-symmetric pertur- 
bations for sufficiently large 7. However, these solutions 
turn out to be always unstable against azimuthal per- 
turbations. In particular, the i/'i" soliton has instabilities 
associated with 1 < m < 5, and the m — 4 mode has the 
largest growth rate for all e. The corresponding eigen- 
value A is real and the eigenfunctions u(r) and v(r) have 
a single maximum near the position of the lateral min- 
imum of the function i?i(r). Following Ref. [jl2| where 
a similar scenario was described for nodal waveguides in 
a saturable self-focusing medium, the above observation 
suggests that the ip^ soliton will decay into 5 solitons 4>o '■ 
one at the origin and four others placed symmetrically 
around it. Next, the t/j^ solution has azimuthal instabili- 
ties with 1 < to < 10. The analysis of the corresponding 
eigenfunctions suggests that, depending on h and 7, the 
decay products will comprise 11 to 13 ^-solitons: 1 at 
the origin; 3 or 4 placed symmetrically around it; and 7 
or 8 forming an outer ring. We verified these predictions 
via direct numerical simulations of the time-dependent 
array (|l|); the simulations corroborated the above sce- 
nario. Thus the nodal solutions can be interpreted as de- 
generate coaxial complexes of the nodeless solitons and 
serve as nuclei of symmetric multisoliton patterns. 

Lastly, we need to understand the stabilisation mech- 
anism in qualitative terms. To this end, we use the vari- 
ational approach. The equation (Q) is derivable from the 
stationary action principle with the Lagrangian 



L = e 27 *Re / {i^* - \Vip\ 2 - |i/>| 2 + \ip\ 4 - h^ 2 )dyi. 
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Choosing the ansatz ip = ^/A e -^~(B+ta)r 2 jTJJfJ with 
A, B, 9, a functions of t, this reduces, in 2D, to 



L B 



h cos( 



0-1 + 

29) cos 



a 

2B 



2B 



A 

cos 2 4> 2 
tan 4> = cr/B. 



(24) 



The 4-dimensional dynamical system defined by (|24|) , has 
two stationary points representing the ip solitons. In 
agreement with the stability properties of the solitons 
in the full PDE, the tp + stationary point is unstable for 
small 7 but stabilises for larger dampings (FigJ^). When 
7 is large we can expand A = Aq + ^Ai + ..., B = Bq+^Bi 
+...,8= 1 + 10! + ..., a= + Letting h = 7 + ^ 
where < c < 1, defining T = - and matching coeffi- 
cients of like powers of yields a 2-dimensional system 



dAo/dT = A [c + 8(Ji 
dBo/dT = 8cti Bo 



H<Ji/B ) 2 } 



2-Bo — jA , 



f 4cti6»i +4((T 2 /B ) 

(7i = ^AqBq 



2B 2 . 



(25) 
(26) 
(27) 



Like their parent system (24), Eqs 



fixed points, the saddle at B Q 
and a stable focus at B, 



(i 



- (J27|) have two 
fc, Aq = ABn 



According to (|5|), the soliton's phase x = f + crr ' con- 
trols the creation and annihilation of the soliton's ele- 
mentary constituents (whose density is \ip\ 2 ). (If Eq.(|^) 
is used as a model equation for Faraday resonance in 
granular media or fluids, J \ip\ 2 dx has the meaning of 
the number of grains or mass of the fluid captured in 
the oscillon.) Since the creation and annihilation occurs 
mainly in the core of the soliton, the variable phase com- 
ponent or 2 plays a marginal role in this process. Instead, 
the significance of the quantity a is in that it controls the 
flux of the constituents between the core and the periph- 
ery of the soliton — see the Xr-term in the r.h.s. of (|). 

If we perturb the stationary point ip + in the 4- 
dimensional phase space of (|24|), the variables 6 and a 
will zap, within a very short time At ~ -, onto the 



2-dimensional subspace defined by the constraints (27). 
After this short transient the evolution of 9 and a will 
be immediately following that of the soliton's amplitude 
\J~A and width 1 / \J~B~. In the case of the ip + soliton, this 
provides a negative feedback: perturbations in A and B 
produce only such changes in the phase and flux that 
the new values of 9 and a stimulate the recovery of the 
stationary values of A and B. (The phase 9 works to re- 
store the number of constituents while a rearranges them 
within the soliton.) In the case of the ip~ the feedback 
is positive: the perturbation-induced phase and flux ( p7j ) 
strive to amplify the perturbation of the soliton's am- 
plitude and width still further. Finally, for small 7 the 
coupling of 9 and a to Aq and -Bo is via differential rather 
than algebraic equations. In this case the dynamics of the 



phase and flux is inertial and their changes may not catch 
up with those of the amplitude and width. The feedback 
loop is destroyed and the soliton destabilises. 
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